% 目标函数
clear all
f = @(x) x(1)^2 + x(2)^2 + x(3)^2 +8;
rho = 0.8;

% 约束条件
% A = [1, -2]
% b = 1
% Aeq = [2,1]
% beq = 1
A = [];
b = [];
Aeq = [];
beq = [];
lb = [];
ub = [];

    %% 这里是碳排的优化_(ADMM形式)
    % 将优化问题通过ADMM来进行拆解（此处针对的碳排最优函数）
    % 构造惩罚函数  名为penalty_function
    % 新目标函数则为 E_all + penalty_function
    f_p_a = @(gama_m) penalty_function(gama_m) + E_all(gama_m);
    gama =  zeros(1, 3);
    gama_1 = zeros(1, 3);
    y = randn(1,3);
    % 采用scaled form形式来进行迭代
    % 首先是u的定义，约束是gama = gama_1
    u=y/rho;
    %然后是二次罚项和拉格朗日乘子的合并
    %alm_item = @(x, z, u) (rho/2) * norm((x - z + u), 2) ^ 2;
    %此处x就是gama，z就是gama_1
    for iteration_tick = 1 : 1000
        k = iteration_tick;
        disp(iteration_tick)

        % 对gama进行优化
        alm_item_gama = @(gama) (rho/2) * norm((gama - gama_1 + u), 2) ^ 2;
        new_gama = @(x) f(x) + alm_item_gama(x);
        [gama,~] = fmincon(new_gama,gama,A,b,Aeq,beq,lb,ub);


        zold = gama_1;
        % 对gama_1进行优化
        alm_item_gama_1 = @(gama_1) (rho/2) * norm((gama - gama_1 + u), 2) ^ 2;

        new_gama_1 = @(x) alm_item_gama_1(x); %penalty_function(x) + alm_item_gama_1(x);
        nonlcon = @opt_power_1;
        gama_1 = fmincon(new_gama_1,gama,A,b,Aeq,beq,lb,ub,'mycon');

        % 对对偶变量进行梯度上升法优化
        u = u + gama - gama_1;
        disp(gama)

    % ----------------------------------------------------------------------判断停止准则---------------------------------------------------------------------------
%          % diagnostics, reporting, termination checks
% 
%         K = gama;    Z = gama_1;
%         
% 
%         history.r_norm(k)  = norm(K - Z);
%         history.s_norm(k)  = norm(-rho*(Z - zold));
% 
%         history.eps_pri(k) = ABSTOL+RELTOL*max(norm(K), norm(-Z));
%         history.eps_dual(k)= ABSTOL+RELTOL*norm(rho*u);
%     %     history.eps_pri(k) = ABSTOL+RELTOL*max(norm(K), norm(-Z));
%     %     history.eps_dual(k)= ABSTOL+RELTOL*norm(rho*u);
% 
% 
%         if ~QUIET
%             fprintf('%3d\t%10.4f\t%10.4f\t%10.4f\t%10.4f\t%10.2f\n', k, ...
%                 history.r_norm(k), history.eps_pri(k), ...
%                 history.s_norm(k), history.eps_dual(k), history.objval(k));
%         end
% 
%     %     if (history.r_norm(k) < history.eps_pri(k) && ...
%     %        history.s_norm(k) < history.eps_dual(k))
%     %          break;
%     %     end

    end





